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ABSTRACT 

The reionization of the Universe, it is believed, occurred by the growth of ionized regions 
(bubbles) in the neutral intergalactic medium (IGM). We study the possibility of detecting 
these bubbles in radio-interferometric observations of redshifted neutral hydrogen (HI) 21 cm 
radiation. The signal (< 1 mJy) will be buried in noise and foregrounds, the latter being at 
least a few orders of magnitude stronger than the signal. We develop a visibility based formal- 
ism that uses a filter to optimally combine the entire signal from a bubble while minimizing 
the noise and foreground contributions. This formalism makes definite predictions on the abil- 
ity to detect an ionized bubble or conclusively rule out its presence in a radio-interferometric 
observation. We make predictions for the currently functioning GMRT and a forthcoming in- 
strument, the MWA at a frequency of 150 MHz (corresponding to a redshift of 8.5). For both 
instruments, we show that a 3 a detection will be possible for a bubble of comoving radius 
Rb > 40 Mpc (assuming it to be spherical) in 100 hrs of observation and Rf, > 22 Mpc in 
1000 hrs of observation, provided the bubble is at the center of the field of view. In both these 
cases the filter effectively removes the expected foreground contribution so that it is below 
the signal, and the system noise is the deciding criteria. We find that there is a fundamental 
limitation on the smallest bubble that can be detected arising from the statistical fluctuations 
in the HI distribution. Assuming that the HI traces the dark matter we find that it will not be 
possible to detect bubbles with Ri, < 8 Mpc using the GMRT and Rj, < 16 Mpc using the 
MWA, however large be the integration time. 
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1 INTRODUCTION 

Quasar absorption spectra (Becker et al., 2001; Fan et al., 2002) 
and CMBR observations (Spergel et al., 2006; Page et al., 2006) 
together imply that reionization occurred over an extended period 
spanning the redshift range 6 < z < 15 (for reviews see Fan, 
Carilli & Keating 2006; Choudhury & Ferrara 2006a). It is cur- 
rently believed that ionized bubbles produced by the first luminous 
objects grow and finally overlap to completely reionize the uni- 
verse (Barkana & Loeb, 2001; Furlanetto, Zaldarriaga & Hernquist, 
2004). In this paper we consider the possibility of detecting these 
ionized bubbles in redshifted 21 cm neutral hydrogen (HI) maps. 

An ionized bubble embedded in neutral hydrogen will appear 
as a decrement in the background redshifted 21 cm radiation. This 
decrement will typically span across several pixels and frequency 
channels in redshifted 21 cm maps. Detecting this is a big chal- 
lenge because the HI signal (~ 1 mJy or lower ) will be buried in 
foregrounds (Shaver et al., 1999; DiMatteo. et al., 2002; Oh, 1999; 
Cooray & Furlanetto, 2004; Santos, Cooray & Knox, 2005) which 



are expected to be at least 2 — 3 orders of magnitude larger. An ob- 
jective detection criteria which optimally combines the entire sig- 
nal in the bubble while minimizing contributions from foregrounds, 
system noise and other such sources is needed to search for ionized 
bubbles. The noise in different pixels of maps obtained from radio- 
interferometric observations is correlated (eg. Thompson, Moran 
& Swenson (1986)), and it is most convenient to deal with visi- 
bilities instead. These are the primary quantities that are measured 
in radio-interferometry. In this paper we develop a visibility based 
formalism to detect an ionized bubble or conclusively rule it out in 
radio-interferometric observations of HI at high redshifts. 

The paper is motivated by the fact that the Giant Metre-Wave 
Radio Telescope (GMRT 1 ; Swamp et al. 1991) which is currently 
functional has a band centered around 150 MHz which corresponds 
to HI at z ~ 8.5. There are several low frequency radio telescopes 
which are expected to become functional in the future (eg. MWA 2 , 
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LOFAR 3 , 21 CMA 4 and SKA 5 ) all of which are being designed 
to be sensitive to the epoch of reionization HI signal. In this paper 
we apply our formalism for detecting ionized bubbles to make pre- 
dictions for the GMRT and for one of the forthcoming instruments, 
namely the MWA. For both telescopes we investigate the feasibility 
of detecting the bubbles, and in situations where a detection is fea- 
sible we predict the required observation time. For both telescopes 
we make predictions for observations only at a single frequency 
(150MHz), the aim here being to demonstrate the utility of our for- 
malism and not present an exhaustive analysis of the feasibility of 
detecting ionized bubbles in different scenarios and circumstances. 
For the GMRT we have used the telescope parameters from their 
website, while for the MWA we use the telescope parameters from 
Bowman et al. (2006). It may be noted that MWA is expected to 
be gradually expanded in phases, and we have used the parame- 
ters corresponding to an early stage, the MWA - Low Frequency 
Demonstrator. 

It is expected that detection of individual bubbles would com- 
plement the studies of reionization through the global statistical sig- 
nal of the redshifted 21 cm radiation which has been studied exten- 
sively (eg. Zaldarriaga, Furlanetto & Hernquist 2004; Morales & 
Hewitt 2004; Bharadwaj & Ali 2005; Bharadwaj & Pandey 2005; 
for a recent review see Furlanetto et al. 2006). 

The outline of the paper is as follows: In Section 2 we dis- 
cuss various sources which are expected to contribute in low fre- 
quency radio-interferometric observation, this includes the signal 
expected from an ionized bubble. In Section 3 we present the for- 
malism for detecting an ionized bubble, and in Section 4 we present 
the results and discuss its implications. The cosmological parame- 
ters used throughout this paper are those determined as the best-fit 
values by WMAP 3-year data release, i.e., Q. m — 0.23, Q. b h 2 = 
0.022, n 3 = 0.96, h = 0.74, a 8 = 0.76 (Spergel et al. 2006). 



2 DIFFERENT SOURCES THAT CONTRIBUTE TO LOW 
FREQUENCY RADIO OBSERVATIONS 

The quantity measured in radio-interferometric observations is the 
visibility V(U, v) which is measured in a number of frequency 
channels v across a frequency bandwidth B for every pair of an- 
tennas in the array. For an antenna pair, it is convenient to use 
U — d/\ to quantify the antenna separation d projected in the 
plane perpendicular to the line of sight in units of the observing 
wavelength A. We refer to U as a baseline. The visibility is related 
to the specific intensity pattern on the sky l v iff) as 

V(U, v) = J d 2 9A{9)I v (9)e 2 " lS0 (1) 

where 9 is a two dimensional vector in the plane of the sky with ori- 
gin at the center of the field of view, and A(9) is the beam pattern 
of the individual antenna. For the GMRT this can be well approxi- 
mated by Gaussian A{9) = e~ e l e ° where 9o ~ 0.6 #fwhm and 
we use the values 2.28° for 9 Q at 153 MHz for the GMRT. Each 
MWA antenna element consists of 16 crossed dipoles distributed 
uniformly in a square shaped tile, and this is stationary with respect 
to the earth. The MWA beam pattern is quite complicated, and it 
depends on the pointing angle relative to the zenith (Bowman et 



3 http://www.lofar.org/ 

4 http://web.phys.cmu.edu/~past/ 
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al., 2007). Our analysis largely deals with the beam pattern within 
1° of the pointing angle where it is reasonable to approximate the 
beam as being circularly symmetric (Figures 3 and 5 of Bowman 
et al. 2007 ). We approximate the MWA antenna beam pattern as 
a Gaussian with 9o = 18° at 153 MHz. Note that the MWA pri- 
mary beam pattern is better modeled as A(ff) oc cos 2 (K 6), but 
a Gaussian gives a reasonable approximation in the center of the 
beam which is the region of interest here. Equation (1) is valid only 
under the assumption that the field of view is small so that it can be 
well approximated by a plane, or under the unlikely circumstances 
that all the antennas are coplanar. 

The visibility recorded in 150 MHz radio-interferometric ob- 
servations is a combination of three separate contributions 



V(U, v) = S(U, v) + N(U, v) + F(U, v) 



(2) 



where S(U, v) is the HI signal that we are interested in, N(U, v) 
is the system noise which is inherent to the measurement and 
F(U , v) is the contribution from other astrophysical sources re- 
ferred to as the foregrounds. Man-made radio frequency interfer- 
ence (RFI) from cell phones and other communication devices are 
also expected to contribute to the measured visibilities. Given the 
lack of a detailed model for the RFI contribution, and anticipating 
that it may be possible to remove it before the analysis, we do not 
take it into account here. 



2.1 The HI signal from ionized bubbles 

According to models of reionization by UV sources, the early 
stages of reionization are characterized by ionized HII regions 
around individual source (QSOs or galaxies). As a first approxima- 
tion, we consider these regions as ionized spherical bubbles char- 
acterized by three parameters, namely, its comoving radius Rb, the 
redshift of its center z c and the position of the center determined by 
the two-dimensional vector in the sky-plane 9 C . The bubble is as- 
sumed to be embedded in an uniform intergalactic medium (IGM) 
with a neutral hydrogen fraction xui- We use r v to denote the co- 
moving distance to the redshift where the HI emission, received at 
a frequency v = 1420 MHz/(l + z), originated, and define r' v = 
d r v /d v. The planar section through the bubble at a comoving dis- 
tance rv is a disk of comoving radius R v = Rb y/l — (Av/Av b ) 2 
where Av = v c — v is the distance from the the bubble cen- 
ter v c in frequency space with v c = 1420 MHz/ (1 + z c ) and 
Avb = Rb /r'„ c is the bubble size in the frequency space. The bub- 
ble, obviously, extends from v c — Avb to v c + Avb in frequency and 
in each frequency channel within this frequency range the image of 
the ionized bubble is a circular disk of angular radius 9 V = R v /r v \ 
the bubble is not seen in HI beyond this frequency range. Under 
such assumptions, the specific intensity of the redshifted HI emis- 



I v {9) = I v xm 



eli-- 



ei 



http://www. skatelescope.org/ 



Av b 

(3) 

where I v = 2.5 x 10 2 ^ (^) ™ the radiation 

background from the uniform HI distribution and Q(x) is the Heav- 
iside step function. 

The soft X-ray emission from the quasar responsible for the 
ionized region is expected to heat the neutral IGM in a shell around 
the ionized bubble. The HI emission from this shell is expected to 
be somewhat higher than I v (Wyithe & Loeb, 2004). We do not 
expect this to make a very big contribution, and we do not consider 
this here. 
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u Av MHz 

Figure 1. Signal from a spherical ionized bubble of comoving radius Figure % signa[ from a spherical ionized bubble of comoving radius 

10 Mpc as a function of baseline U for different frequency channels. 1Q Mpc as a f^^on of Au = u _ Vc for different baselines. 



If we assume that the angular extent of the ionized bubble is 
small compared to the angular scale of primary beam ie. 9 V <C 80, 
we can take A(9) outside the integral in eq. (1) and write the signal 
as A(9 C ) J d 2 9I„(9)e 27Tieu , which essentially involves a Fourier 
transform of the circular aperture fl— [ 9 — C \ r v /R v ^. For 
example, a bubble of radius as large as 40 Mpc at z — 8.5 would 

o 

have an angular size of only 9 V ~ 0.25 which satisfies the condi- 
tion 9 V <C 9q. In a situation where the bubble is at the center of the 
field of view, the visibility is found to be 

"2Ji(27r<70, 



Scenter(£/, v ) = —TtI v Xm9 v 



2ttU9 v 



9 1- 



Av b 
(4) 

where Ji(x) is the first order Bessel function. Note that 
5ccntcr (U, v) is real and it is the Fourier transform of a cir- 
cular aperture. The uniform HI background also contributes 
I v n9 e to the visibility, but this has been dropped as it is 
quite insignificant at the baselines of interest. Note that the approxi- 
mations used in eqs. (4) have been tested extensively by comparing 
the values with the numerical evaluation of the integral in eq. (1). 
We find that the two match to a high level of accuracy for the situ- 
ations of interest here. In the general situation where the bubble is 
shifted by 9 C from the center of the field of view, the visibility is 
given by 



S(U,u) 



S centei {U,v) (5) 
- e c/ e o drop in the 



i.e., there is a phase shift of e 2 ' !rlU ' <: am j a e _ 
overall amplitude. 

Figures 1 and 2 show the U and Av dependence of the vis- 
ibility signal from an ionized bubble with Rb = 10 Mpc located 
at the center of the field of view at v c — 153 MHz (z c = 8.3), 
assuming xm = 1. The signal extends over Av = ±Avi where 
Avb — 0.56 MHz. The extent in frequency Avb = Rb/r' Vc scales 
oc Rb when the bubble size is varied. The Bessel function Ji(x) 
has the first zero crossing at x — 3.83. As a result, the sig- 
nal S(U,v) extends to U = 0.61r„[R by /l - {Av/Avb) 2 ]' 1 
where it has the first zero crossing, and Uq scales with the bub- 
ble size as Uo oc 1/ Rb- The peak value of the signal is S(0, v) = 
■Ki v (Rbjr v f sj\ - (Av/Av b ) 2 and scales as S(Q,u) oc R 2 b if 
the bubble size is varied. We see that the peak value of the signal 
is S(0, v c ) = 70 /iJy for bubble size Rb = 10 Mpc and would 
increase to 1.75 mjy if Rb = 50 Mpc. Detecting these ionized 
bubbles will be a big challenge because the signal is buried in noise 
and foregrounds which are both considerably larger in amplitude. 
Whether we are able to detect the ionized bubbles or not depends 



critically on our ability to construct optimal filters which discrimi- 
nate the signal from other contributions. 



2.2 HI fluctuations 

In the previous sub-section, we assumed the ionized bubble to be 
embedded in a perfectly uniform IGM. In reality, however, there 
would be fluctuations in the HI distribution in the IGM which, in 
turn, would contribute to the visibilities. This contribution to the HI 
signal can be treated as a random variable S(U, v) with zero mean 
{S(U, v)) — 0, whose statistical properties are characterized by 
the two- visibility correlation (S(Ui,vi)S(U2,V2))- This is related 
to Pni(k) the power spectrum of the 21 cm radiation efficiency in 
redshift space (Bharadwaj & Ali, 2004) through 



(S(U 1 ,v)S*(U 2 ,v + A„))=5 fIiA ^ 

r 00 

x / dfc|| .Fkt(k) co&{k\\r' v Au) 
Jo 



(6) 



where 5ft n is the Kronecker delta ie. different baselines are un- 
correlated, To estimate the contribution from the HI fluctuations 
we make the simplifying assumption that the HI traces the dark 
matter, which gives Pm(k) = x\u (1 +/i 2 ) P(k) where P(k) 
is the dark matter power spectrum and fi is the cosine of the angle 
between k and the line of sight. This assumption is reasonable be- 
cause the scales of interest are much larger than the Jeans length 
A j ~ 10 — 100 kpc, and we expect the HI to cluster in the same 
way as the dark matter. 

In addition to the above, there could be other contributions to 
the HI signal too. For example, there would be several other ionized 
regions in the field of view other than the bubble under consider- 
ation. The Poisson noise from these ionized patches will increase 
the HI fluctuations and there will also be an overall drop in the 
contribution because of the reduced neutral fraction. These effects 
will depend on the reionization model, and the simple assumptions 
made in this paper would only provide a representative estimate 
of the actual contribution. Figure 3 shows the expected contribu- 
tion from the HI fluctuations (HF) to the individual visibilities for 
GMRT and MWA. Note that while this can be considerably larger 
than the signal that we are trying to detect (particularly when the 
bubble size is small), there is a big difference between the two. 
The signal from the bubble is correlated across different baselines 
and frequency channels whereas the contribution from random HI 
fluctuations is uncorrelated at different baselines and it become un- 
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correlated beyond a certain frequency separation Av (Bharadwaj 
& Ali, 2005; Datta, Choudhury & Bharadwaj, 2006). 



2.3 Noise and foregrounds 

The system noise contribution N(U, v) in each baseline and fre- 
quency channel is expected to be an independent Gaussian random 
variable with zero mean ({N) = 0) and whose variance is indepen- 
dent of U and v c . The predicted rms. noise contribution is (Thomp- 
son, Moran & Swenson (1986)) 



2005; Datta, Choudhury & Bharadwaj, 2006) 



(N*) 



A eff y/Av c At 



(V) 



where T sys is the total system temperature, fcs is the Boltzmann 
constant, A e f j is the effective collecting area of each antenna, Av c 
is the channel width and At is correlator integration time. Equation 
(7) can be rewritten as 



(iV 2 ) = C x 



Av c 



_y 1/a (*LY 1/a 

1MHz ) V lsec J 



(8) 



where C x varies for different interferometric arrays. Using the 
GMRT parameters T sys = 482K and A eff /2k B = 0.33K/Jy 
at 153MHz gives C x = 1.03Jy for the GMRT where as for MWA 
T S ys = 470K and A eff /2k B = 5 x 10~ 3 K/Jy (Bowman et al., 
2006) gives C x = 65.52Jy. The rms noise is reduced by a factor 
^At/tots if we average over t obs /At independent observations 
where t a b s is the total observation time. Figure 3 shows the ex- 
pected noise for a single baseline at 153MHz for Av c — 50 KHz 
and an observation time of 100 hrs for both the GMRT and MWA. 
Though Tsys is nearly equal for the GMRT and the MWA, the noise 
in a single baseline is expected to be 60 times larger for MWA 
than that for the GMRT. This is a because the individual antennas 
have a much larger collecting area at the GMRT as compared to the 
MWA. The fact that the MWA has many more antennas (N = 500) 
as compared to the GMRT (N = 30) compensates for this. Note 
that nearly half (16) of the GMRT antennas are at very large base- 
lines which are not particularly sensitive to the signal on the an- 
gular scales the ionized bubble, and only the other 14 antennas in 
the 1 km x 1 km central square will contribute towards detecting 
the signal. For both the GMRT and the MWA, T sys is dominated 
by the sky contribution T s k y with the major contribution coming 
from our Galaxy. We expect T 3ys to vary depending on whether 
the source is in the Galactic plane or away from it. The value which 
we have used is typical for directions off the Galactic plane. Fur- 
ther, the noise contribution will also be baseline dependent which 
is not included in our analysis. 

Contributions from astrophysical foregrounds are expected to 
be several order of magnitude stronger than the HI signal. Extra- 
galactic point sources and synchrotron radiation from our Galaxy 
are predicted to be the most dominant foreground components. As- 
suming that the foregrounds are randomly distributed, with possible 
clustering, we have {F(U, v)) = for all the baselines other than 
the one at zero spacing (U = 0), which is not considered in this pa- 
per. The statistical properties are characterized by the two- visibility 
correlation {F(Ui, vi)F(U2, ^2)). We express this (details in Ap- 
pendix A) in terms of the multi-frequency angular power spectrum 
(hereafter MAPS) Ci(fi V2) of the brightness temperature fluctu- 
ations at the frequencies v\ and V2 as (Santos, Cooray & Knox, 



'dB\ (dB\ 

G2irU 1 (Vl V2)- 



dT \dT, 



(9) 



where {dB/dT) v = 2k B v 2 /c 2 is the conversion factor to spe- 
cific intensity, and we have assumed that the primary beam pattern 
A{9) — e~ B is frequency dependent through 0o oc z/ _1 and 
use 0i and 02 to denote the value of 0q at v\ and 1/2 respectively. 
Note that the foreground contribution to different baselines are ex- 
pected to be uncorrelated. 

For each component of the foreground the MAPS is modeled 



(10) 



where Vf = 130 MHz, and for each foreground component A, 
(3 and a are the amplitude, the power law index of the angular 
power spectrum and the mean spectral index respectively. The ac- 
tual spectral index varies with line of sight across the sky and this 
causes the foreground contribution to decorrelate with increasing 
frequency separation Av = \v\ — vi\ which is quantified through 
the foreground frequency decorrelation function V2) (Zaldar- 
riaga, Furlanetto & Hernquist, 2004) which has been modeled as 



h(y\ V2) = exp 



log? (^)/2C 2 



(11) 



We consider the two most dominant foreground components 
namely extragalactic point sources and the diffuse synchrotron ra- 
diation from our own galaxy. Point sources above a flux level 
Scut can be identified in high-resolution continuum images and 
removed. We note that absence of large baselines at the MWA re- 
stricts the angular resolution, but it may be possible to use the large 
frequency bandwidth 32 MHz to identify continuum point sources 
in the frequency domain. S C ut depends on a the rms. noise in the 
image. We use S C ut = 5a where a is the rms noise in the image 
given by (assuming 2 polarizations) 

c ( b \- 1/2 /u 5 y 1/2 

^/2Nt V 1MHz J V l sec / 



where Nb = N(N — l)/2 is the number of independent baselines, 
N is the number of antennas in the array, B is the total frequency 
bandwidth and t b s the total observation time. For t bs = 100 hrs 
and B = 6 MHz we have Scut = 0. lmjy for the GMRT and using 
B = 32 MHz it gives S cut = 0.2 mjy for the MWA. The value of 
Scut will be smaller for longer observations, but reducing S cu t any 
further does not make any difference to our results so we hold S cu t 
fixed at these values for the rest of our analysis. The confusion noise 
from the unresolved point sources is a combination of two parts, the 
Poisson contribution due to the discrete nature of these sources and 
the clustering contribution. The amplitude of these two contribu- 
tions have different S C ut dependence. The parameter values that 
we have used are listed in Table 1 . We have adopted the parameter 
values from Santos, Cooray & Knox (2005) and incorporated the 
Scut dependence from DiMatteo. et al. (2002). 



Figure 3 shows the expected foreground contributions for the 
GMRT and MWA. The galactic synchrotron radiation is the most 
dominant foreground component at large angular scales (U < 1000 
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Figure 3. The magnitude of the different contributions to the visibility V(U, v) at v = 153 MHz as a function of U. The signal, foregrounds (FG), noise (NS) 
and HI fluctuations (HF) contributions are shown for the GMRT (left) and MWA (right). The expected signal is shown for bubbles with radius R = 10 Mpc 
and R = 50 Mpc. The noise is estimated for a single baseline assuming an observation time t (, s = 100 hrs and channel width Av c = 50 KHz. 



Table 1. Parameters values used for characterizing different foreground 
contributions 



Foregrounds 




j4(mK 2 ) 


a 




5 


Galactic synchrotron 




700 


2.80 


2.4 


4 


Point source 


61 


( s„ t \ - 5 


2.07 


1.1 


2 


VO.lmJyJ 


(clustered part) 












Point source 


0.16 


V0.1mJ y y 


2.07 





1 


(Poisson part) 













for GMRT and U < 2000 for MWA), while the clustering of the 
unresolved extragalactic point sources dominates at small angular 
scales. For all values of U, the foregrounds are at least four orders 
of magnitude larger than the signal, and also considerably larger 
than the noise. 

The MWA has been designed with the detection of the statis- 
tical HI fluctuation signal in mind, and hence it is planned to have 
a very large field of view. The foreground contribution to a single 
baseline is expected to be 10 times stronger for the MWA than 
for the GMRT because of a larger field of view. As we shall show 
later, the increased foreground contribution is not a limitation for 
detecting HII bubbles. The foregrounds have a continuum spectra, 
and the contribution at two different frequencies at a separation Av 
are expected to be highly correlated. For Av = 1 MHz, the fore- 
ground decorrelation function Ii(Av) falls by only 2 x 10~ 6 for 
the galactic synchrotron radiation and by 3 x 10~ 5 for the point 
sources. In contrast, the signal from an ionized bubble peaks at a 
frequency corresponding to the bubble center and falls rapidly with 
Av (Figure 2). This holds the promise of allowing the signal to be 
separated from the foregrounds. 



3 FORMALISM FOR DETECTING THE IONIZED 
BUBBLE 

We consider a radio-interferometric observation of duration t b s , 
carried out over the frequency range v\ to v-z. The HI signal from 
an ionized bubble, if it is present in the data, will be buried in fore- 
grounds and noise both of which are expected to be much larger. In 
this Section we present a filtering technique aimed at detecting the 



signal from an ionized bubble if it is present in our observations. To 
detect the signal from an ionized bubble of radius Rt with center at 
redshift z c (or frequency v c ) and at an angle 9 C from the center of 
the field of view, we introduce an estimator E[Rt, z c , 9 C ] defined 

as 



E : 



(13) 



where Sf(U,u) is a filter which has been constructed to detect 
the particular ionized bubble. Here U a and vt refer to the different 
baselines and frequency channels in our observations, and in eq. 
(13) we are to sum over all independent data points (visibilities). 
Note that the estimator E and the filter Sf(U, v) both depend on 
[Rb, z c , 9 C ], the parameters of the bubble we wish to detect, but we 
do not show this explicitly. The values of these parameters will be 
clear from the context. 

We shall be working in the continuum limit where the two 
sums in eq. (13) can be replaced by integrals and we have 



E = 



d 2 U 



dvp N {U,v) Sf*(U,v)V{U,v) (14) 



d 2 Udv pn(U,v) is the fraction of data points ie. baselines and 
frequency channels in the interval d 2 U dv. Note that pn{U,v) 
is usually frequency dependent, and it is normalized so that 
J d 2 U J dv Pn(U, v) = 1. We refer to Pn(U, v) as the normal- 
ized baseline distribution function. 

We now calculate (E) the expectation value of the estimator. 
Here the angular brackets denote an average with respect different 
realizations of the HI fluctuations, noise and foregrounds, all of 
which have been assumed to be random variables with zero mean. 
This gives (V0, v)) = S(U, v) and 



(E) 



dvp N (U,v) S f *(U,v)S(U,v) 



(15) 



We next calculate the variance of the estimator which is the 
sum of the contributions from the noise (NS), the foregrounds(FG) 
and the HI fluctuations (HF) 

{(AE) 2 ) ee {(E~{E)) 2 ) 

- ((A^) ns + ((A^) fg+ ((A^) hf . 

(16) 

To calculate the noise contribution we go back to eq. (13) and use 
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the fact that the noise in different baselines and frequency channels 
are uncorrelated. We have 



((AE) 2 ) NS = (TV ) 



E I S f&> 



''hi 



(17) 



which in the continuum limit is 



((A£) 2 ) NS 



(TV 2 )/ El 



d 2 U dvp N (U,v) \S f (U,v)\ 2 (18) 



The term y |^(TV 2 ) / ^2 a b lj is the same as a, the rms. noise in the 
image (eq. 12). We then have 

((A£) 2 ) NS =a 2 jd 2 U f dup N {U,v) \S f (U,u)\ 2 . 

(19) 



For the foreground contribution we have 

(AS) 2 



fg = J d 2 [/i J d 2 U 2 J dv! J du 2 

X p N (Ul, Vl)PN{Un,V2)Sf*(Ui, Vl)Sf*(U 2 , V 2 ) 

X {F(U 1 ,u 1 )F(U 2 , v-z)) (20) 
In the continuum limit we have (details given in Appendix A) 

X C 2itUi {V\,V2) 

which gives the variance of the foreground contribution to be 

dB \ f dB 
dT .. \dT 



I-, 

(21) 



dui / du 2 

\ it i i , . . , 

x pjv((7, vi)p N (U, v 2 )S f *(U, vi)Sf(U, va) 
xC 2 rfh,i/2) (22) 



We use eq. (22) to calculate ( (A_E) 2 ^ 



/ 1 



too, with the differ- 



ence that we use the power spectrum C 2v viy, v + Av) for the HI 
fluctuation from Datta, Choudhury & Bharadwaj (2006) instead of 
the foreground contribution. 

In an observation it will be possible to detect the presence of 
an ionized bubble having parameters [Rb,Zc,0c] at, say 3-sigma 



confidence level, if (E) > 3y {(AE) 2 ). In such a situation, an 
observed value E can be interpreted as a detection with 99.7% 



(i.e., 3-sigma) confidence if E > 3y {(AE) 2 ). The presence of 
the ionized bubble can be ruled out at the same level of confidence 



if (E)-E > 3J((AE) 2 ) 



3.1 Baseline distribution 

In this subsection we discuss the normalized baseline distribution 
function pn(U,v) which has been introduced earlier. Figure 4 
shows the baseline coverage for 14hrs of observation towards a 
region at declination 8 = 45° with the GMRT at 153MHz. In this 
figure u and v refer to the Cartesian components of the baselines 
U. Note that the baseline distribution is not exactly circularly sym- 
metric. This asymmetry depends on the source declination which 
would be different for every observation. We make the simplifying 
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Figure 4. This shows the baseline coverage for 14 hrs of GMRT 153 MHz 
observation at 45° declination. 




Figure 5. This shows the normalized baseline distribution pn(U,v) for 
the GMRT and the MWA at 153 MHz. The wiggly curve shows the actual 
values for the GMRT observation shown in Figure 4 and the smooth curve 
is the analytic fit. 



assumption that the baseline distribution is circularly symmetric 
whereby Pn{U, v) is a function of U. This considerably simplifies 
our analysis and gives reasonable estimates of what we would ex- 
pect over a range of declinations. Figure 5 shows pn{U , v) for the 
GM RT determined from the baseline coverage shown in Figure 4. 
We find that this is well described by the sum of a Gaussian and 
an exponential distribution. The GMRT has a hybrid antenna dis- 
tribution (Chengalur et al., 2003) with 14 antennas being randomly 
distributed in a central square approximately 1 km x 1 km and 16 
antennas being distributed along a Y each of whose arms is 14 km 
long. The Gaussian gives a good fit at small baselines in the central 
square and the exponential fits the large baselines. Determining the 
best fit parameters using a least square gives 



p N {U,v) 



1 

B 



A 



lkm 



0.21 exp 



+ 



9.70 x 10" 3 exp 



U\ 
2a 2 
UX-b 



d 



(23) 



where a = 0.382 km, b = 0.986 km, d = 3.07 km and B is the 
frequency bandwidth which has a maximum value of 6 MHz. 

Following Bowman et al. (2006) we assume that the MWA an- 
tennas are distributed within a radius of 0.750 km with the density 



of antennas decreasing with radius r as p a nt(r) 



and with a 
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maximum density of one antenna per 18m 2 . The normalized base- 
line distribution is estimated in terms of p ant (r) and we have 



Pn{U,v) = 



4.4 x 

p a „t(\r- \U\)d(p 



t(r) 



(24) 



the 

„2 



bandwidth 

r2 \2 



B 



32 MHz, 



where 

\r-\U\ = (r 2 + U 2 A 2 - 2r A U cos<f) l/2 . Note that p Ar ({7, v) 
depends on the observed frequency. Figure 5 shows the normal- 
ized baseline distribution function pn(U,u) for both the GMRT 
and the MWA. We see that maximum baseline for the GMRT is 
Umax ~ 10, 000 whereas Umax ~ 750 for the MWA. However, 
the smaller baselines will be sampled more densely in the MWA as 
compared to the GMRT. 



3.2 Filter 

It is a major challenge to detect the signal which is expected to be 
buried in noise and foregrounds both of which are much stronger 
(Figure 3). It would be relatively simple to detect the signal in a 
situation where there is only noise and no foregrounds. The signal 
to noise ratio (SNR) is maximum if we use the signal that we wish 
to detect as the filter (ie. S f (U, v) = S(U, v)) and the SNR has a 
value 



((E)) 



1 



NS 



dvp N (U,v) | S(U,u) 



(25) 



The observing time necessary for a 3-cr detection (i.e., SNR = 3) 
would be the least for this filter. The difficulty with using this fil- 
ter is that the foreground contribution to \J ((AE) 2 ) is orders of 

magnitude more than ((E)). The foregrounds, unlike the HI sig- 
nal, are all expected to have a smooth frequency dependence and 
one requires filters which incorporate this fact so as to reduce the 
foreground contribution. We consider two different filters which 
reduce the foreground contribution, but it occurs at the expense of 
reducing the SNR, and t obs would be more than that predicted by 
eq. (25). 

The first filter (Filter I) subtracts out any frequency indepen- 
dent component from the frequency range v c — B /2 to v c + B /2 
with B' < B ie. 

S f (U, v )=^j a [s(U, v) 

B Jur-B'/2 



(26) 



This filter has the advantage that it does not require any prior 
knowledge about the foregrounds except that they have a contin- 
uous spectrum. It has the drawback that there will be contributions 
from the residual foregrounds as all the foregrounds are expected 
to have a power law spectral dependence and not a constant. A 
larger value of B causes the SNR to increases, and in the limit 
B — > oo the SNR approaches the value given in eq. (25). Unfor- 
tunately the residues in the foregrounds also increase with B . We 
use B — 4Ai/f, provided it is less than B, and B = B otherwise. 



The frequency dependence of the total foreground contribu- 
tion can be expanded in Taylor series. Retaining terms only up to 
the first order we have 

Ci(u u i/ 2 ) = Ci(v c , v c ) [1 - (Ai/i + Ai/ 2 ) a e ff/v c ] (27) 



where At 



v c and a 



eff 



E,a. (iooo/o^ 18 the ettec " 



tive spectral index, here i refers to the different foreground com- 
ponents. Note that a e / / is I dependent. The second filter that we 
consider (Filter II) allows for a linear frequency dependence of the 
foregrounds and we have 



S f (U,v) = (1 + Q e //Ai//z/ c ) 



A c 



[s(U, v) 



6(1-2 | v-v c | IB') p+ s '/2 



Jv c -B' /2 



S(U,v')dv 



(28) 



Note that for both the filters we include an extra factor (A c /A) . 
This is introduced with the purpose of canceling out the A 2 depen- 
dence of of the normalized baseline distribution function Pn(U, v) 
and this substantially reduces the foreground contribution. 



4 RESULTS AND DISCUSSIONS 

We first consider the most optimistic situation where the bubble 
is at the center of the field of view and the filter center is exactly 
matched with the bubble center. The size distribution of HII re- 
gions are quite uncertain, and would depend on the reionization 
history and the distribution of ionizing sources. However, there 
are some indications in the literature on what could be the typi- 
cal size of HII regions. For example, Wyithe et al. (2005) deduce 
from proximity zone effects that Rb « 35 Mpc at z £s 6, which 
should be considered as a lower limit. On the other hand, Furlan- 
etto et al. (2006) (Figure 1(a)) infer that the characteristic bubble 
size R b > 10 Mpc at z = 8 if the ionized fraction Xi > 0.75 
(Rb ~ 50Mpc if Xi ~ 0.9). Theoretical models which match a 
variety of observations (Choudhury & Ferrara, 2007) imply that Xi 
could be as high as 90 % at 2 ~ 8, which would mean bubble 
sizes of ~ 40 — 50 Mpc. To allow for the large variety of possibil- 
ities, we have presented results for a wide range of Rb values from 
2 Mpc to 50 Mpc. We restrict our analysis to a situation where the 
IGM outside the bubble is completely neutral (a'Hi = 1). The sig- 
nal would fall proportional to xm if the IGM outside the bubble 
were partially ionized (xhi < 1). The expected signal (E) and 3- 

sigma fluctuation 3 x <J {(AE) 2 ) from each of the different com- 
ponents discussed in Sections 2 and 3 as a function of bubble size 
Rb are shown in Figures 6 and 7. Both the figures show exactly the 
same quantities, the only difference being that they refer to Filter I 
and Filter II respectively. A detection is possible only in situations 

where (E) > 3 x \J ((AE) 2 ), the rhs. now refers to the total 
contribution to the estimator variance from all the components. 

The signal is expected to scale as i?f and the noise as R^ 2 in 
a situation where the baseline distribution is uniform ie. Pn(U, v) 
is independent of U. This holds at U < 300 for the GMRT (Figure 
5), and the expected scaling is seen for R b > 20 Mpc. For smaller 
bubbles the signal extends to larger baselines where Pn(U, v) falls 
sharply, and the signal and the noise both have a steeper R b depen- 
dence. The MWA baseline distribution is flat for only a small U 
range (Figure 5) beyond which it drops. In this case the signal and 
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R b Mpc R b Mpc 

Figure 6. The signal quantified through the expectation value of the estimator (E) for Filter I. The other components (NS - Noise, FG - Foregrounds, HF - HI 
Fluctuations) are quantified through their contribution to the 3-sigma fluctuation 3 X y ((A_E) 2 ). 
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Figure 7. Same as Figure 6 except that Filter II is used instead of Filter I. 
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Figure 8. The observing time t j, s that would be required for a 3 a detection 
of a bubble of radius R b provided it is at the center of the field of view. 
The vertical lines shows the lower limit (due to HI fluctuations) where a 
detection will be possible (R b = 8 Mpc for GMRT and Rb = 16 Mpc for 
MWA). 

noise are found to scale as and Rl respectively. Note that the 
maximum baseline at MWA is U = 750, and hence a considerable 
amount of the signal is lost for Rb < 10 Mpc. 

At both the GMRT and the MWA, for 100 hrs of observation, 
the noise is larger than the signal for bubble size Rt < 40 Mpc. At 
the other extreme, for an integration time of 10 5 hrs the noise is be- 
low the signal for R b > 6 Mpc for the GMRT and R b > 8 Mpc for 



the MWA. The foreground contribution turns out to be smaller than 
the signal for the entire range of bubble sizes that we have consid- 
ered, thus justifying our choice of filters. Note that Filter II is more 
efficient in foreground subtraction, but it requires prior knowledge 
about the frequency dependence. For both the filters the foreground 
removal is more effective at the GMRT than the MWA because of 
the frequency dependence of Pn(U, v). The assumption that this 
is proportional to A 2 is valid only when p(U, v) is independent of 
U, which, as we have discussed, is true for a large U range at the 
GMRT. The A dependence is much more complicated at the MWA, 
but we have not considered such details here as the foreground con- 
tribution is anyway smaller than the signal. It should also be noted 
that the foreground contribution increases at small baselines (eq. 
10), and is very sensitive to the smallest value of U which we set 
at U = 20 for our calculations. Here it must be noted that our re- 
sults are valid only under the assumption that the foregrounds have 
a smooth frequency dependence. A slight deviation from this and 
the signal will be swamped by the foregrounds. Also note that this 
filtering method is effective only for the detection of the bubbles 
and not for the statistical HI fluctuations signal. 

The contribution from the HI fluctuations impose a lower limit 
on the size of the bubble which can be detected. However long be 
the observing time, it will not be possible to detect bubbles of size 
Rb < 8 Mpc using the GMRT and size Rb < 16 Mpc using the 
MWA. The HI fluctuation contribution increases at small baselines. 
The problem is particularly severe at MWA because of the dense 
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A 6 arc min A 8 arc min 

Figure 9. The Overlap between the signal and the filter when there is a mismatch A8 between the centers of the bubble and the filter for GMRT (left) and 
MWA (right). The results are shown for different bubble sizes. 



sampling of the small baselines and the very large field of view. We 
note that the MWA is being designed with the detection of the sta- 
tistical HI fluctuation signal in mind, and hence it is not surprising 
that this contribution is quite large. For both telescopes it may be 
possible to reduce this component by cutting off the filter at small 
baselines. We have not explored this possibility in this work be- 
cause the enormous observing times required to detect such small 
bubbles makes it unfeasible with the GMRT or MWA. 

Figure 8 shows the observation time that would be required to 
detect bubbles of different sizes using Filter I for GMRT and the 
MWA. Note that the observing time shown here refers to a 3 a de- 
tection which is possibly adequate for targeted searches centered on 
observed quasar position. A more stringent detection criteria at the 
5 a level would be apropriate for a blind search. The observing time 
would go up by a factor of 3 for a 5 a detection. The observing time 
is similar for Filter II and hence we do not show this separately. In 
calculating the observing time we have only taken into account the 
noise contribution as the other contributions do not change with 
time. The value of Rb below which a detection is not possible due 
to the HI fluctuations is shown by vertical lines for both telescopes. 
We see that with 100 hrs of observation both the telescopes will 
be able to detect bubbles with Rb > 40 Mpc while bubbles with 
Rb > 22Mpc can be detected with 1000 hrs of observation. 

The possibility of detecting a bubble is less when the bubble 
centre does not coincide with the centre of the field of view. In fact, 
the SNR falls as e _e< =^ e ° if the bubble center is shifted away by 9 C 
from the center of the field of view and the filter is also shifted so 
that its center coincides with that of the bubble. There will be a cor- 

ofl 2 / o 2 

responding increase t a b s oc e <= ' in the observing time required 
to detect the bubble. It will be possible to detect bubbles only if they 
are located near the center of the field of view (9 C <C 9o), and the 
required observing time increases rapidly with 9 C for off-centered 
bubbles. 

When searching for bubbles in a particular observation it will 
be necessary to consider filters corresponding to all possible value 
of Rb, v c and 9 C . A possible strategy would be to search at a dis- 
crete set of values in the range of Rb, v c and 9 C values where a 
detection is feasible. The crucial issue here would be the choice 
of the sampling density so that we do not miss out an ionized 
bubble whose parameters do not exactly coincide with any of the 
values in the discrete set and lie somewhere in between. To illus- 
trate this we discuss the considerations for choosing and optimal 
value of A9 C the sampling interval for 9 C . We use {E)[A9] to de- 
note the expectation value of the estimator when there is a mis- 
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Figure 10. Same as the Figure 8 considering three different antenna distri- 
butions pant(r) oc 1/r, 1/r 2 , 1/r 3 for the MWA. 

match A9 between the centers of the bubble and the filter. The ra- 
tio Overlap = (E) [A9]/(E)[0], shown in Figure 9 for GMRT 
(left panel) and MWA (right panel), quantifies the overlap between 
the signal and the filter as A9 is varied. We see that the choice of 
A9 would depend on the size of the bubble we are trying to detect 
and it would be smaller for the GMRT as compared to the MWA. 
Permitting the Overlap to drop to 0.9 at the middle of the sampling 
interval, we find that it is 8' at the GMRT and 20' at the MWA for 
R b = 50 Mpc. 

The MWA is yet to be constructed, and it may be possible 
that an antennae distribution different from p an t( r ) oc 1/r 2 , may 
improve the prospects of detecting HII bubbles. We have tried out 
Pant(r) oc 1/r and 1/r 3 for which the results are shown in Fig- 
ure 10. We find that the required integration time falls consider- 
ably for the 1/r 3 distribution whereas the opposite occurs for 1/r. 
For example, for Rb = 50 Mpc the integration time increases by 
5 times for 1/r and decreases by 3 times for 1/r 3 as compared 
to 1/r 2 . Based on this we expect the integration time to come 
down if the antenna distribution is made steeper, but this occurs at 
the expense of increasing the HI fluctuations and the foregrounds. 
We note that for the 1/r 3 distribution the foreground contribution 
is more than the signal, but it may be possible to overcome this 
by modifying the filter. The increase in the HI fluctuations is in- 
evitable, and it restricts the smallest bubble that can be detected to 
R = 26MPc for 1/r 3 . In summary, the 1/r distribution appears 
to be a good compromise between reducing the integration time 
and increasing the HI fluctuations and foregrounds. 
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Finally we examine some of the assumptions made in this 
work. First, the Fourier relation between the specific intensity and 
the visibilities (eq. 1) will be valid only near the center of the field 
of view and full three dimensional wide-field imaging is needed 
away from the center. As the feasibility of detecting a bubble away 
from the center falls rapidly, we do not expect the wide-field effects 
to be very important. Further, these effects are most significant at 
large baselines whereas most of the signal from ionized bubbles is 
in the small baselines. 

Inhomogeneities in the IGM will affect the propagation of ion- 
ization fronts, and the ionized bubbles are not expected to be ex- 
actly spherical (Wyithe et al., 2005). This will cause a mismatch be- 
tween the signal and the filter which in turn will degrade the SNR. 
In addition to this, in future we plan to address a variety of other 
issues like considering different observing frequencies and making 
predictions for the other upcoming telescopes. 

Terrestrial signals from television, FM radio, satellites, mobile 
communication etc., collectively referred to as RFI, fall in the same 
frequency band as the redshifted 21cm signal from the reionization 
epoch. These are expected to be much stronger than the expected 
21cm signal, and it is necessary to quantify and characterize the 
RFI. Recently Bowman et al. (2007) have characterized the RFI for 
the MWA site on the frequency range 80 to 300MHz. They find 
an excellent RFI environment except for a few channels which are 
dominated by satellite communication signal. The impact of RFI 
on detecting ionized bubbles is an important issue which we plan 
to address in future. 

The effect of polarization leakage is another issue we postpone 
for future work. This could cause polarization structures on the sky 
to appear as frequency dependent ripples in the foregrounds inten- 
sity . This could be particularly severe for the MWA. 
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APPENDIX A: RELATION BETWEEN 

VISIBILITY- VISIBILITY CORRELATION AND MAPS 

In this appendix we give the calculations for expressing the two vis- 
ibility correlation in terms of the Multi-frequency angular power 
spectrum (MAPS). We can write the visibility V(U, v) as a 
two-dimensional Fourier transform of the brightness temperature 
T(8, v) [see equation (1)] 

V( P'" )= {m) J d2eA (0,^T(e,u)e 2 ^° (Al) 

where (8B/dT) v is the conversion factor from temperature to spe- 
cific intensity and A(6, v) is the beam pattern of the individual an- 
tenna. The visibility- visibility correlation is then given by 

(V { U^)V02^))=(^) (If 



x j d 2 e J d 2 e' 'A{e,u{)A{e > 2 ) 

x (T(9, fi)T(0\ V2 ^ e ^<S-u l+ S'-u 2 ) 



(A2) 

The correlation function for the temperature fluctuations on the 
sky would simply be the two-dimensional Fourier transform of the 
MAPS C2ttU {vi , v?) 



{T{0, vi)T{9' , V2)) = J d 2 UC 2 „u(vi,V2)t 



-2m(§-6')-U 



(A3) 



Using the above equation in equation in (A2), we obtain 



x j d 2 U C 2 7vu{vi,v 2 ) 

x A{U 1 -U,v 1 )A02 + U,v2) (A4) 
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where A(U, v) is the Fourier transform of the beam pattern 
A(Q, v). If the beam pattern is assumed to be Gaussian A(8, v) = 
e -° / e o > the Fourier transform too is given by a Gaussian function 

A(U,v) = n9? ) e- 7T2u20 " (A5) 
Hence, the visibility correlation becomes 

'dB\ fair 

. dT J { dT . 



{V{Ui,ui)V(Ui,V2)) 



I 



7v 2 [(u 1 -u) 2 ej +(u 2 +&) 2 e%] 



(A6) 



where 8i and 82 are the values of #0 at v\ and V2 respectively. 
Now, since the two Gaussian functions in the above equation is 
peaked around different values of U, the integrand will have a non- 
zero contribution only when \Ui + fa | < (ir max[9i, fo]) -1 . 
In case the typical baselines are much larger than the quantity 
(n max[0i, fe]) -1 , the integral above can be well approximated 
as being non-zero only when Ui = — U2. Then 



<V(C/i,^i)\/((7 2 ,^)) 



dB 
dT 



dB 
dT 



X TV 6 1 8 2 C2ttU 1 {v\,Vl) 



/■ 



2 [(U 1 -U) 2 (e{ + 2 )] 



8\8\ 



rt+ei 



dB 
dT 



dB_ 

dT 



X C 2 nU 1 (l'l, V2) 



"2 

(A7) 



which is what has been used in equation (9). 

In the continuum limit, the Gaussian A(U, v) can be approx- 
imated by a delta function, i.e., A(U, v) w 5$ (U) (which corre- 
sponds to the limit 8 a — > 00); the visibility-visibility correlation is 
then given as 

<V(tfi,«*)V(tf 2) i*)> = 5^0, + U2) (||) 

X.C2i,U 1 {vi,V2) (A8) 



which corresponds to equation (21) in the main text. 



